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ABSTRACT 

We present a characterization of the local luminosity-color, bi-variate distribution of IRAS galaxies from the 
1.2 Jy sample, selected at 60 /im. The TZ(60, 100) infrared color is used as the best single parameter description 
of the IR spectral energy distribution of galaxies. We derive an analytical form of the distribution and use it to 
constrain the effect of the IR color distribution on evolution models for high redshift, farlR luminous galaxies. 
Our adopted evolution retains the locally observed correlation between luminosity and color, such that the larger 
characteristic luminosities at higher redshift have a warmer characteristic color The width of the color distribution 
at a given luminosity remains constant for all redshifts. We demonstrate that there is the potential for both hotter 
and colder sources to be missed in cosmological surveys. An evolving bi-variate luminosity function coupled 
with the cold source bias of sub-mm selected surveys suggests the existence of a large population of cold sources 
appearing in such surveys. Likewise, a hot source bias for most SIRTF wavebands together with a bi-variate model 
suggests an excess of hot sources being selected. We test the evolutionary form against available data for higher 
redshift, farlR galaxies. The data do not reveal evidence for any strong evolution in the characteristic luminosity- 
color distribution as a function of redshift over < z < 1 . However, there is marginal evidence for a broadening 
of the color distribution at higher redshifts, consistent with our locally characterized trend of a broadening in the 
IR color distribution at the highest luminosities. 



Subject headings: galaxies: evolution — galaxies: formation 
galaxies 



sub-mm: galaxies — radio: galaxies — IR: 



1. INTRODUCTION 

The local infrared-luminous galaxies detected by the IRAS 
satellite exhibit a vast array of source properties. Infrared color 
has typically been used to attempt to parametrize the IRAS pop- 
ulation (Soifer & Neugebauer 1991). Dale et al. (2001) demon- 
strated that the Seo^m/Sioo/jm fux ratio or color (hereafter re- 
ferred to as 7?.(60, 100)) provides the best single parameter char- 
acterization of the IRAS galaxy population, in addition to lumi- 
nosity. While a complicated array of dust properties contribute 
to the spectral energy distribution (SED) of each galaxy, stud- 
ies of IRAS galaxies have typically reduced the description to a 
best fit single dust temperature, Td, with a one-to-one mapping 
to 7?.(60, 100). Indeed, changing the dust temperature has been 
demonstrated to have a significantly larger effect on the galaxy 
SED than dust emissivity, mid-IR spectral index, or cosmol- 
ogy (Blain et al. 2002). The inferred luminosity of an infrared 
galaxy for a fixed observed mid-IR flux density increases by a 
factor of 10 if the dust temperature is doubled. 

It has been demonstrated that low-redshift IRAS galaxies 
exhibit slowly varying correlations of 71(60, 100) with lumi- 
nosity (Dale et al. 2001; Dunne et al. 2000; Andreani & 
Franceschini 1996). However, over a large spread of lumi- 
nosities the 7^^(60, 100) ratio of infrared galaxies does change 
systematically. Fitting single dust temperature (Td) models to 
7e(60,100), we find -20 K for low-redshift spirals (Reach et 
al. 1995; Alton et al. 2000; Dunne & Eales 2001) and 30 K- 
60 K for the high luminosity objects typically detected by IRAS 
(Soifer & Neugebauer 1991, Stanford et al. 2000). High red- 
shift, hyper-luminous galaxies can show dust temperatures of 
up to 1 10 K (e.g., Lewis et al. 1998), implying a continuation in 
the luminosity-Td relation out to the higher luminosities char- 



acteristic of the distant Universe. 

However, while a statistical relation exists between the 
7^(60, 100) and IR luminosity, the distribution is broad. 
We find in substantial numbers both extremely luminous, 
yet cold galaxies, as well as low luminosity, hot galax- 
ies. One surprisingly cold and luminous galaxy, Arp 302 
(or UGC9618/NGC5051), has been identified from the IRAS 
bright galaxy sample (BGS) with S6o^m=6.8 Jy, Si()0pni=15.3 Jy, 
and Lfir=3.89x 10" L0. This is the system with largest devia- 
tion from the median TZ{60, 100) for its luminosity. Lo, Gao & 
Guendl (1997) suggest this galaxy is the most massive known 
(in terms of CO gas mass). By contrast, galaxies with very 
hot IR color, and without obvious AGN contributions, have 
been identified from the faintest Lfir sources detected by IRAS. 
NGC1377, NGC4491 and IRAS1953 with 7e(60,100)-l and 
Lfir~ 10^ L0 are as hot or hotter than the Ultra-Luminous In- 
frared Galaxy (ULIG), Arp220, but with 10"^ the FIR luminos- 
ity (Roussel et al. 2003). 

The possible importance of cold, luminous galaxies to farlR 
and sub-mm surveys has been pointed out by Eales et al. (1999, 
2000). This has recently been highlighted observationally by 
Chapman et al. (2002a) who demonstrated that cold and lu- 
minous sources exist at higher redshift, identifying two ISO- 
FIRBACK sources (FBI -40 and FBI -64) with Ltir> lO^^ L0 
with zll galaxies, and finding best fit single temperature grey- 
bodies of 26 K and 31 K respectively (dust emissivity f3 = 1.6). 
The local cold source, Arp 302, approximately matches the 
color and luminosity of FBI -40. However, FBI -40 and FBl- 
64 were discovered from a random sampling of the sub-mm 
luminous FIRBACK sources, and they cannot be considered as 
an insignificant portion of the high-z ULIG population. 



2 



The bi-variate distribution of IRAS galaxies 



It is therefore a concern that studies of the evolving galaxy 
populations typically assume a small range of template galaxy 
spectral energy distributions (SEDs), sampling only the mono- 
tonic relation of dust temperature to luminosity in choosing the 
SED templates (Blain et al. 1999a,b, Malkan & Stecker 2000, 
Rowan-Robinson 2001, Chary & Elbaz 2001, Chapman et al. 
2002b, Franceschini et al. 2002). Modeling efforts have in ef- 
fect tied the dust temperature directly to the farlR luminosity, 
ignoring the Td distribution for each luminosity class. This in 
part may be a result of the fact that the TZ(6Q, 100) distribution 
for local IRAS galaxies has never been carefully studied, and no 
analytical form for the distribution has been presented. 

Importantly, an evolving distribution bi-variate in luminos- 
ity and color, C), consistent with the broad distribution 
observed locally should subsume a non-negligible fraction of 
both cold, luminous galaxies and hot, faint galaxies. Surveys 
which select objects at either the cold Raleigh- Jeans tail of the 
dust SED, or the hot Wien tail, will preferentially detect ap- 
propriately cold or hot objects for a given luminosity class, if 
they exist in non-negligible numbers. For instance, the cold 
source bias of farlR (e.g., 170/im) and submillimeter selected 
surveys will result in any existing population of cold sources 
being over-represented. 

At the highest redshifts, the poorly understood submillime- 
ter (sub-mm) population is thought to dominate the most lumi- 
nous IR galaxies. Since the discovery of high redshift sub-mm 
sources (Small, Ivison & Blain 1997), there has been an ongo- 
ing debate about their nature and their dust properties. In the 
absence of redshifts, it is difficult to understand whether they 
represent hot sources at very high redshifts, or colder sources 
at more modest redshifts (e.g. Eales et al. 1999,2000), similar 
to FBI -40 and FBI -64 (Chapman et al. 2002a). Chapman et al. 
(2003a) have measured spectroscopic redshifts for radio iden- 
tified sub-mm galaxies, claiming a typical dust temperature of 
~40 K by assuming the empirical relation between the far-IR 
and radio observed locally (e.g., Helou et al., 1985). However, 
without additional SED measurements, it is unclear what distri- 
bution in dust properties for the sub-mm galaxies remains con- 
sistent with their currently measured properties. 

In this paper, we characterize the local luminosity-color, bi- 
variate distribution, $(£, C), of IRAS galaxies from the 1.2Jy 
sample, selected at 60 /im. We represent the infrared luminos- 
ity with Ltir, and the 7<!.(60, 100) IR-color is used as the best 
single parameter description of the IR SED of galaxies. Ltir is 
defined as in Dale et al. (2001), integrating over the SED from 
3-1 100 /im. These authors define a bolometric conversion be- 
tween the more typical FIR luminosity (e.g., Helou et al., 1985) 
as 

log(TIR / FIR) = ao + aix + a2X^ + aT,x^ + a4X* 

where x = log ^jj^^y and [a(z=0] = [0.2738, -0.0282, 0.7281, 
0.6208, 0.9118]. The first part of the paper is directed at the 
derivation of an analytical form for the local IRAS color dis- 
tribution. In order to constrain the effect of the IR color dis- 
tribution on our understanding of high redshift farlR luminous 
galaxies, we evolve the bi-variate distribution according to lu- 
minosity evolution prescriptions used in the literature. For con- 
venience, we use the terms 'color' (meaning 60/100 micron flux 
ratio) and 'temperature' interchangeably, even though we do 
not collapse the SED description to a single temperature black- 
body. The color is taken to be a diagnostic of the typical heat- 
ing conditions in the ISM of a galaxy, and therefore would be 
indicative of a characteristic dust temperature. Section 2 de- 
scribes our sample used to construct $(£, C). Section 3 pro- 



vides a statistical analysis of the 1 .2 Jy sample, and provides a 
best fit analytical form to $(£, C). Section 4 explores the evo- 
lutionary behavior of $(£, C), while Section 5 compares the 
evolving <&(£, C) to existing data sets at higher redshift. 

2. SAMPLE SELECTION 

Our starting point for this study is the Seo^m > 12 Jy sam- 
ple of galaxies (Fisher et al. 1995). We plot the 7^(60, 100) 
distribution of these galaxies as a function of both Lfir and 
Ltir in Fig. 1. Lfir is calculated directly from the 60 /im and 
100/im flux densities and the redshift as described in Helou et 
al. (1985). The Ltir parameter represents a bolometric correc- 
tion for the flux from 3 to 1 100 microns, and represents a larger 
correction for sources with cooler TZ(6Q, 100) colors, as seen 
in the comparison of the two panels of Fig. 1. The analytical 
expression mapping the average FIR luminosity to TIR lumi- 
nosity from Dale et al. (2001) was reproduced in the previous 
section. 

While studying this sample, we noticed a significant tail of 
very cold and luminous galaxies, departing bi-modally from the 
main distribution (see Fig. 1). Upon closer inspection of the 60 
and 100/im IRAS fluxes using the XSCANPI application pro- 
vided by the Infrared Processing and Analysis Center (IPAC), 
we found that almost 100% of these sources have spurious mea- 
surements. Firstly, the cirrus contamination can be large, affect- 
ing the 100/im more than 60 /im and leading to an apparently 
small 60/im/lOO/im ratio. Secondly, the 60 /im flux measure- 
ments found for some of the sources from XSCANPI were less 
than 1 .2 Jy, and should not be included in our catalog. 




,*» ■ • — : 

■ « 




ll»-<> 



Fig. I. — The Lfir and Ltir distributions in 7^.(60, 100). Spurious cold, 
luminous sources liave been removed from the sample, as indicated here by 
solid lines delimiting the region of spurious sources. The medians in equal 
number bins are plotted as large circles, with the interquartile range shown in 
error bars. 

We flagged all such sources in our sample and removed them 
from subsequent consideration. To ensure that this effect did 
not significantly contaminate the main distribution, we ran- 
domly selected 200 of the sources with Lfir > 10" for XS- 
CANPI analysis, finding them to have correct flux estimates 
to within 5% of the Fisher et al. (1995) values. We then derive 
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the TIR luminosity function with the revised catalog of 1 .2Jy 
sources. We adopt an accessible volume technique as described 




log (Ltir/Lq) 
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with $(L)AL as the number density of sources (Mpc"-') in 
the luminosity range L to L+ AL. The accessible volume, V,, 
represents the /th source in the sample, the maximum volume 
in which the object could be located and still be detected in the 
IRAS 1.2Jy 60 /im catalog. The sum is then over all sources 
within the luminosity range. We then map sources to their TIR 
luminosity using the above definition. Our constructed TIR lu- 
minosity function is therefore directly related to the 60 fim lu- 
minosity function for this sample. 

Figure 1 also shows the median 7^(60, 100) values and the in- 
terquartile range, revealing the dust temperature to luminosity 
relation. In the subsequent sections, we shall study and charac- 
terize this relation. 
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3. RESULTS 
3.1. Statistics 

In choosing a luminosity variable in C) to describe the 
distribution of galaxies, it is desirable to minimize the depen- 
dencies between C and C. Ltir has been shown to provide less 
dependency than Lfir on IRAS colors (Dale et al. 2001). Simi- 
larly, 7?,(60, 100) has been demonstrated to parametrize the vari- 
ation in mid-farlR properties better than any other combination 
of IRAS bands (Dale et al. 2001). 



Fig. 2. — Upper Panel: The 1st quartile, median, and last quartile for the 
1.2 Jy sample with a log scaling in 7?.(60, 100), suggesting a roughly constant 
width and shape in log space. Lower Panel: The same statistics with the me- 
dian subtracted, shown in a linear scaling emphasizing the small Hnear increase 
in width with increasing luminosity. The width is constant in log(7?.(60, 100)). 

Having chosen our variable representation of $(£, C), we 
can proceed to quantify the relation for the S6o>1.2Jy IRAS 
sample and explore the implications. We begin by calculating 
the median and first/last quartile statistics for the 1.2Jy sam- 
ple, corrected for spurious cold, luminous sources as described 
in section 2. The analytic relation in the distribution median is 
well fit by a dual power law, as shown in figure 2: 

log7^(60/100) = 0.162 xlog(LTiR/i0)-2.O8O, Lur >5 x 10'°Lc 

log 7e(60/100) = 0.022 x log(LTiR/i0)- 0.593, Ltir < 5 x 10'"L, 
We express this composite function with a smooth transition 
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as: 



7^(60/100) = C* x(l + — ^r''x(l + -^)'' 



with,7 = 0.16, 6 = 0.02, C* = 0.45, = 5.0 x 10^% 



The statistics and the fitted relation are shown in Fig. 2a. We 
also fit the width of the inter-quartile statistics after removing 

the median (Fig. 2b). The width of the distribution described in 
terms of this statistic is essentially constant in log(7<l(60, 100)) 
as a function of Ltir, therefore showing a slight broadening 
in 7?.(60, 100) as shown in Fig. 2b. At the highest luminosi- 
ties there is some evidence for significant broadening of even 
the log(7?.(60, 100)), although the statistics are poor and there 
remains a concern that all spurious sources have not been re- 
moved from the sample (section 2). 

These equations provide a simple, first order description of 
the $(£, C) distribution of IRAS galaxies from which to gage a 
more comprehensive analysis. This expression should be suit- 
able for many apphcations and evolutionary extrapolations. 




0.25 0.5 0.75 

^60 (U.m/^100 iJ.m 

Fig. 3. — The C) distribution in log IR-color, sliown as a smoothed 
histogram. The effective histogram bin size is d7?.(60, 100)=0.05. Histograms 
of LxiR are shown offset by a fixed constant for clarity, and labeled on the right 
in L0. 



3.2. Fitting the ^(C, C) distribution 

We now consider the detailed distributions in 7?,(60, 100) for 
the 1.2 Jy sample. Histograms of the distribution from each lu- 
minosity bin are plotted in Fig. 3. For clarity, a fixed constant is 
introduced to offset each class of Ljir. With no apriori assump- 
tions about the 72.(60, 100) distribution of the IRAS sources, we 
begin by testing whether the population might be well repre- 
sented by a gaussian in the 7^(60, 100) color. A gaussian in 
both Unear Tl{60, 100) and log7^(60, 100) was first fit to each 
bin of 380 sources. The functional form 



C+B*x+D* exp[-l /2 x (^-i^f] 



where x represents 7?.(60, 100), was fit to each histogram of 380 
sources. An equal weighting of the points in the histogram was 
assumed in the fitting. The residuals were compared, with the 
being a factor ^^2 smaller for the gaussian in log (or log 
normal distribution), suggesting the log normal as the more ap- 
propriate representation. 
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Fig. 4. — upper panel: The TZ(60, 100) distribution (heavy solid lines) with 
log normal fits (light green lines) and linear Gaussian fits (light dashed ma- 
genta lines). Offsets for each luminosity class are inserted as in Fig. 3. The 
median 7?.(60, 100) has been removed for comparison of luminosities, bottom 
panel: The subtracted residuals are shown for both log normal (heavy line) 
and linear gaussian fits (dashed light line). Note that an improved fit, in the 
sense of reduced cannot be obtained by a first order term of the form x*G. 
A higher order term is required to obtain a better fit, which is not justified for 
the present data set. 

The resuhs of our fitting are shown in Fig. 4a, with each lu- 
minosity bin of the color histogram offset for clarity. The me- 
dian has been removed from each histogram in order to com- 
pare directly the form of the distributions at various values of 
Ltir. Fits of log normal distribution (light green line) and linear 
gaussian (light dashed magenta line) are shown. The residuals 
are plotted in Fig. 4b, emphasizing the superior fit of the log 
normal, as well as the skewed nature of the distribution. Heavy 
lines show the residuals from a log gaussian fit, while the light 



dashed magenta lines show the residuals from a linear gaussian 
fit. While the skew is an apparently large and systematic effect, 
(it appears as a double-S, integral shape) the most significant 
deviations occur in the wings of the distribution, where our er- 
ror bars are larger than the skew due to the rapidly diminishing 
number of sources in each histogram bin. For completeness 
we attempted to fit the skewed component with a multiplicative 
factor jc X exp[-z^/2], but that did not improve the x^- The sys- 
tematic skew to the distribution therefore cannot be removed 
without a term which is too high an order to justify from the 
statistics of the sample. We therefore use the log-normal distri- 
bution fit as the best representation of the 1 .2 Jy IRAS galaxy 
sample. Note also that the residual skew is negligible near the 
intermediate luminosity bins of the population (the best sam- 
pled region) and becomes more pronounced towards the lumi- 
nosity extremes. The skew may therefore partially be a result of 
including broader luminosity ranges in the distribution to retain 
equal number bins. 

The distribution over 7^(60, 100) color can then be expressed 

as 



g(C) = exp(-l/2x[(C-C,)/acf), 



where ac = 0.065, C = log7e(60/100), 



C, = C.x(l + -^)-^x(l + ^)^ 



with, 7 = 0.16, 5 = 0.02, C* = 0.45, = 5.0 x W% 



As the residual skew is close to symmetrical about the gaus- 
sian fit peak, we find that our fit to the distribution describes 
nearly the same function as the median and first/last quartile ex- 
pression (a log-normal) presented in the previous section. The 
formal $(£, C) is then expressed as follows, where we have de- 
rived the TIR luminosity function directly from the 1 .2Jy cata- 
log as described in the previous section, 

<^>(C,C) dCdC = $!(£) X $2(C) dCdC 



= P*x4- ^ X (1 + ^r'' X exp[-l/2(^)2]d/:dC, 



with£ = LTiR, p*=4.34M;9c-^L-\ = lO'^'^L©, a=1.55, /3 = 2.10, 



C = log(56o^,„/5ioo,™), C, as above, crc =0.065. 
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Fig. 5. — The range of local SEDs used in our model from the catalog of 
Dale et al. (2001, 2002). The flux units are ergs/s/cm^/Hz normalized to the 
median luminosity of a = 0. 1 galaxy in our model. 

Our analysis has suggested that the C-C distribution is best 
described with a dual power-law in C with a break luminosity 
5.0 X 1O'"L0, a faint end slope of J = 0.02 and a bright end 
slope of 7 = 0.16. We offer two possible physical interpreta- 
tions of this result. 

Firstly, this may be the point of transition from cirrus- 
dominated luminosity to active star formation (i.e. high-density 
photo-dissociation regions, etc.) dominating the luminosity. 
The different power laws then arise because in the former 
case the luminosity increases mostly by making the emitting 
dust mass larger, whereas for active star formation, the heating 
drives the luminosity. 

Secondly, this may be the point of transition from where the 
60 /im band still has fluctuating grain emission contributing, to 
having the 60 /^m band dominated by large dust grains. Once 
large dust grains dominate the 6Qfim flux, the 7?,(60, 100) ra- 
tio begins to look like a black-body ratio, leading, naturally to 
the broken power-law relation. However in this case, the high 
luminosity, steep portion of the relation should scale with dust 
temperature like CxT^. Since this is not observed, this second 
explanation cannot dominate the observed relation. 

4. CONSEQUENCES OF THE C-C RELATION FOR EVOLUTION 
MODELS 

Several authors have recently modeled the evolution of dusty 
galaxies using pure, or nearly pure, luminosity evolution, repro- 
ducing the observed counts and backgrounds at IR through sub- 
mm wavelengths (e.g. Blain et al. 1999a,b; Malkan & Stecker 
2000; Rowan-Robinson 2001; Chary & Elbaz 2001; Chapman 
et al. 2002b; Franceschini et al. 2001). None of these models 
include a bi-variate luminosity function (LF), and at best map 
the dust temperature or spectral shape monotonically to the lu- 
minosity. 

The spectral shape of the far-infrared background (FIRB) de- 
tected by DIRBE at 140 and 240 fim (Puget et al. 1996; Fixsen 
et al. 1998) indicates a peak at ^ 200 /im. However, the width 



of the peak suggests that galaxies over a large range in redshifts 
and/or dust temperature contribute to the FIRB. There is a de- 
generacy of dust temperature with redshift (Blain 1999), since 
both translating a source to higher redshift and decreasing the 
dust temperature will shift the SED to lower frequency. This 
suggests an urgency to explore the effect of the bi-variate LF 
on evolutionary models. 

4. 1 . A color added evolution model 

Our goal is to understand the key differences in <!>(£, C) from 
the single variable <&(£). We turn to Monte Carlo simulations 
of a pure luminosity evolution paradigm, similar to the models 
of the above authors. In companion papers (Lewis et al. 2003, 
S. Chapman in preparation), we demonstrate that our model is 
able to simultaneously fit the farlR background and the counts 
at various wavelength bands. 

We evolve the local FIR luminosity function (LF) using 
$(L, v) = ^o(L/g{z), 1^0(1 +z)). Our evolution function follows 
a power law in redshift, g(z) = (1+z)"^, out to z=2.6. Be- 
yond z=2.6 the function drops again as g{z) = (1 +2)""* to avoid 
over-predicting the farlR background. This power-law index 
is chosen based on evolutionary models fit to both optical and 
sub-mm wavelength data (Blain et al. 1999a,b). The evolu- 
tion to z = 2.6 was chosen as this provided the best fit to the 
joint radio/sub-mm sample of galaxies (Chapman, Lewis & 
Helou 2002, Lewis et al. 2003) - currently the best constraint 
on high-z farlR galaxy evolution (Blain et al. 1999a). The red- 
shift z = 2.6 is also close to the consensus of the median red- 
shift for SCUBA sources suggested for example by photomet- 
ric redshifts in Ivison et al. (2002). In order to match the pure 
luminosity evolution, we truncate the evolved functions at low 
luminosity and at both extremes of color, so that they integrate 
to identical numbers of sources per comoving volume. 

One additional assumption in this model is that the bi-variate 
distribution, $(£, C), evolves in £ as a function of redshift, 
while the relation as observed locally between C and C holds 
at all redshifts. The width of the distribution in C follows our 
measured relation from section 3, which remains close to con- 
stant in log(7?,(60, 100)). This implies that the characteristic 
dust temperature of higher redshift sources will be greater be- 
cause of the rise in C 

We note however, that while the no evolution hypothesis for 
C is the simplest, it is not necessarily the obvious choice physi- 
cally. Higher redshift sources may have distinctly different dust 
properties than local analogs. For instance, metallicities will 
likely be lower, there may be less dust to heat, and dust may be 
less centrally concentrated. This form of evolution might lead 
to the TZ{6Q, 100) distribution being biased to cooler IR colors 
with increasing redshift. In the absence of any strong constraint 
from the counts or background (Lewis et al. 2003, S. Chapman 
in preparation), this assumption needs to be tested by future ob- 
servational data sets. However, we shall examine in the follow- 
ing section the preliminary evidence in support of this simple 
hypothesis. 

We then incorporate this hypothesis into our Monte Carlo 
models, drawing 7^(60, 100) values from the evolving $(£, C) 
to fill our desired survey volume. This model does not include 
observational error, but should accurately represent the chance 
of observing a source with a given IR color in a particular sur- 
vey. The simulation is cut off at 10** L© in order to avoid trun- 
cation effects from the fixed comoving number density in the 
luminosity functions. 

In order to tie our evolving £(TIR) luminosity function to ob- 
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Fig. 6. — A plot of the L-Tj plane showing points drawn from our Monte Carlo sampling of the evolving bi-variate LF, and depicting all objects in the simulated 
volume of the survey regardless of their detectability by various instruments. For comparison we show a representation of the same evolving model with a very 
narrow range of colors for each luminosity (bold, offset points). The absissa has been mapped from 7^(60, 100) to dust temperature (1^) using the single temperature 
greybody which provides the best fit to the SED template. The quantized appearance of the points in Tj results from the individual 64 SED classes considered. The 
LF with narrow range of IR color closely approximates a single variable LF whereby luminosity maps monotonically to Tj or IR color 



servable quantities, we map each galaxy with a given Ltir and 
TZ{6Q, 100) value to a spectral energy distribution (SED) shape. 
Normalization of the SEDs is accomplished by integrating the 
SED over the 3-1 100 /im range and scaling to the adopted Ltir. 
Spectral templates are taken from the Dale et al. (2001, 2002) 
catalog, divided into 64 classes from U(60, 100)=0.29 to 1.64, 
corresponding roughly to single component dust temperature 
models of 19 K to 56 K. Representative SEDs from our cata- 
log are shown in Fig. 7. The global galaxy SEDs adopted here 
use a power-law distribution of dust over heating intensity U in 
order to reproduce the range of photometric and spectroscopic 
properties observed by IRAS and ISO for galaxies in the local 
Universe. Dale & Helou (2002) specifically investigate the dust 
emissivity of the local IRAS galaxies, using the sub-mm data 
presented in Dunne et al. (2000). They find the dust emissivity 
should vary with environment with /3 = 2.5-0.41ogZ//, where U 
is the local radiation field. 

4.2. Model Results and Implications 

To underscore the effect of not including the $(£, C) distri- 
bution in the redshift evolution, we show in Figs. 6-8 the out- 



put from our luminosity evolution models, with and without the 
broad local color distribution in the LF. 

In Fig. 6, we plot the C-T^ plane showing points drawn from 
our Monte Carlo sampling of the evolving bi-variate LF. We 
compare directly the bi-variate $(£, C) to an equivalent $(£), 
with a very narrow range of colors for each luminosity (bold, 
offset points). The ordinate has been mapped from 7^(60, 100) 
to dust temperature (Td) using the single temperature greybody 
which provides the best fit to the SED template. The quantized 
appearance of the points in Td results from the individual 64 
template SED classes considered. The LF with narrow range 
of IR color closely approximates a single variable LF whereby 
luminosity maps monotonically to Tj or IR color. Fig. 6 pro- 
vides an overview of 7?,(60, 100) as a function of C, simply ex- 
trapolating the local C-C relation to the required luminosities 
- our baseline assumption in this model. The distribution is 
symmetric in log(Td) about each luminosity. As the $(£, C) 
distribution remains identical to our locally characterized form 
(§ 3) for all redshifts, any deviations from the median local re- 
lation reflect only the scatter in the relatively small numbers of 
luminous galaxies in this simulated survey volume. 
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Fig . 7 . — Compai'ison of our Monte Cai'lo representation of the evolving C) to an evolving model LF with a one-to-one mapping of luminosity to color. The 
lighter symbols represent only luminosity and thus follow identical distributions in both cases. The difference between models is apparent when flux limited surveys 
are defined on the model points. We overlay flux limited samplings of the model points at 24 /im and 850/im with darker symbols. Structure in the 24 /jm survey 
is a result of PAH bands (rest ~l()/^m) being redshifted through the SIRTF 24/jm filter. In the single variable LF model, luminosities map uniquely to fluxes and 
lines can be drawn which reflect fixed temperatures, shown at 50 K, 38 K and 35 K, top to bottom. 



In Fig. 7, we show the Ljir distribution as a function of red- 
shift. The left panels show the $(£, C) model, while the right 
panels show the single variable LF model. Until we apply a 
flux limit to the figure (darker points), there is no difference 
between the visualizations since they have the same luminosity 
evolution formalism. Fig. 7 is our Monte Carlo representation 
of the evolving LF; vertical slices reveal the dual power law 
$(£) at each redshift. Note that with comoving number density, 
in the absence of luminosity evolution, a local physical density 
of N objects per Mpc^ corresponds to a physical density of N * 
(1 at redshift z. In this type of model the evolution is not 



necessarily meant to imply that ULIRGs fade into less lumi- 
nous LIRGs with time. Instead, it provides a picture where, for 
instance, the local counterparts with volume densities similar to 
ULIRGs at z ^ 1 have lO'^'^L© at the present-day. This form 
of evolution is then not describing the internal physics of the IR 
luminous population; high-z ULIRGs do not evolve into similar 
local ULIRGs, but instead likely form ellipticals (e.g., Sanders 
& Mirabel 1996; Tacconi et al. 2002), fading away from the LF. 

When flux limited surveys are considered in the context of 
Fig. 7 (850 /im - top, 24/_tm - bottom), differences in the two 
models become manifest. In the case of the single variable dis- 
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Fig. 8. — Histograms of sources with redsliift are plotted for the cold, medium and hot range in IR color (again mapped to the equivalent single component dust 
temperature). Heavy dashed lines represent the evolving C) distribution, while lighter solid lines represent the single variable LF. An excess of hotter sources 
by a factor >2 is revealed for the ^(C, C) model where all sources are luminous enough to be detected regardless of the observing band. For the medium and cold 
temperature slices, the effective luminosity ranges lie near the sensitivity limits of the survey which additionally affect the histograms. 



tribution, each Ltir point maps uniquely to a flux for a given 
wavelength. However in the bi-variate LF, each Ltir point cor- 
responds to a probability distribution of fluxes corresponding 
to the log-normal distribution in 71(60, 100) and the associated 
range of SED templates that can be tied to the Ltir value. The 
sensitivity limits in the 850 fim and 24 /im bands are shown for 
characteristic survey depths with SCUBA and SIRTF (2mJy 
r.m.s., and 0.1 mJy r.m.s., respectively). The structure in the 
24 /im survey is a result of PAH bands (rest '^lO/im) being 
redshifted through the SIRTF 24 /im filter. At 850 /im, the K- 
correction from the rising grey-body dust SED results in a flat 
Luminosity-Flux relation for redshifts z Zl- 

The effect is subtle in Fig. 7, as both the luminosity func- 
tion and the color distribution are scattering the observed fluxes, 
largely canceling dramatic differences in the effective luminos- 
ity limits probed with redshift. However, the most important 
difference between the bi-variate and luminosity-only models 
is apparent in Fig. 7: in the simpler $(£) model, the flux limit 
for a given wavelength translates at each redshift into a transi- 
tion range of luminosities, within which galaxies are or are not 
detected depending on their color We draw lines on the single 
variable model for C corresponding to Td of 35 K, 38 K, and 



50 K, using the mapping illustrated in Fig. 6. 

The implications of this model difference are therefore much 
more apparent in a direct analysis of Td, as demonstrated in 
Fig. 8. Histograms of sources versus redshift are plotted for the 
850 //m and 24 /im bands, comparing the detectability of hot- 
ter and cooler IR colors (again mapped to the equivalent single 
component dust temperature). Heavy dashed lines represent the 
evolving $(£, C) distribution, while lighter lines represent our 
approximation to a single variable LF. 

We first consider the model differences for a SCUBA 850 fim 
survey in Fig. 8. A cut for detectable sources colder than 
Td <35K (shown in Fig. 7) lies along the sensitivity limit of 
SCUBA. This results in the most dramatic difference between 
the models, since in the bi-variate model, the steep luminosity 
function scatters many more cold, under-luminous sources into 
the sample than the corresponding loss of hotter, higher lumi- 
nosity sources. The selection along the cold side of the grey- 
body peak for SCUBA therefore produces a substantial popu- 
lation of cold sources which are strikingly absent when only 
the single variable LF is included in the evolution. This differ- 
ence results in an important prediction of the existence of cold 
and luminous galaxies at moderate to high redshifts - a predic- 
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tion which has been verified through the detection of ULIRGs 
with Td < 30K (Chapman et al. 2002a). The opposite effect is 
shown in the middle panel of Fig . 8 for the range 35<Tci<38K, 
producing a band in luminosity lying just above the sensitivity 
limit of our SCUBA survey in Fig. 7. The bi-variate model now 
scatters the warmer sources below our detection limit, resulting 
in a dearth of sources relative to the single variable LF. For the 
hottest sources (Tj > 50 K), the luminosities are large enough 
(Fig. 7) that the sensitivity limit of the survey has no bearing 
on the observed sources, and the comparison of models sim- 
ply reflects the steep luminosity function which asymmetrically 
scatters more low luminosity hot sources into the temperature 
cut than luminous cold sources which fall beneath the cut. The 
bi-variate $(£, C) model thus predicts almost twice as many 
warm galaxy detections, even though the underlying luminos- 
ity distribution is the same as in the $(£) model. 

For surveys selecting sources along the hot dust side of the 
grey-body peak, that being the case for all the accessible wave- 
lengths of SIRTF except 170/im, colder luminous sources will 
be missed and hotter low luminosity sources will be preferen- 
tially detected. This is exactly in the opposite sense to the situa- 
tion for sub-mm selected surveys (such as with SCUBA), pref- 
erentially detecting colder sources at a given luminosity (Eales 
et al. 1999, Blain et al. 2002). 

The sensitivity limit of SIRTF is adversely affected by the 
steep Wien slope of the far-mid-infrared SED, making more 
distant sources difficult to detect and resulting in the steep sen- 
sitivity curve in Fig. 7. The result in Fig. 8 is that the histograms 
near the sensitivity limit of SIRTF will be dominated by the 
large numbers of sources lying both above and below any fixed 
temperature cut. The histograms in the cold and warm 24 /^m 
panels therefore continue to show the asymmetric shuffling of 
sources due to the steep luminosity function, but only as a small 
perturbation on the large number of sources which are unaf- 
fected by the sensitivity limit. However, for hotter dust tem- 
peratures, an excess of a factor greater than two of sources is 
predicted by the $(£, C) model. This excess occurs for the 
SIRTF 24fim band for the same reason discussed above for 
SCUBA. All sources are luminous enough to be detected re- 
gardless of their temperature, and there are far more lower lu- 
minosity sources which are boosted by the bi-variate LF into 
the Td > 50 K bin than vice versa. 

Consideration of the model differences in Fig. 8 also ad- 
dresses the question of the SCUBA/SIRTF galaxy overlap. As 
the SCUBA and SIRTF observing windows lie on opposite 
sides of the dust spectral energy peak, the direct overlap of the 
populations is a strong function of the dominant dust tempera- 
ture. In particular, a cold (Tj < 30 K) SCUBA population may 
be difficult to detect even in the deepest SIRTF exposures, as 
revealed in the high redshift excess of cold SCUBA sources 
compared with the SIRTF population. Scrutiny of Fig. 7 along 
the survey flux hmits reveals that at z * 3, the individual SCUBA 
galaxies detected in either model are mostly detected by SIRTF 
at 24 /xm at this survey depth. However, deeper sub-mm surveys 
by instruments like the Large Millimeter Telescope (Hughes et 
al. 2001) may uncover fainter cold galaxies at lower redshifts 
which would require very deep SIRTF exposures to detect. By 
contrast, most or all of the SCUBA sources should be detected 
by SIRTF over all redshifts when the characteristic dust temper- 
atures are higher, as suggested by the similar N(z) in the lower 
panels of Fig. 8. 

Figs. 7&8 are also suggestive of the effect of the evolutionary 
form on the detectabiUty of galaxies at different wavelenghs. 



The properties of the sources detectable at high redshift are a 
sensitive function of the evolutionary form adopted - in our 
scenario, lower redshift peaks in the evolution function would 
result in an appreciably smaller fraction of high redshift galax- 
ies detected (effectively generating differences in the redshift 
distributions). For example, a direct application of our $(£, 
C) model to the SCUBA population (Chapman, Lewis & Helou 
2002, Lewis et al. 2003) has suggested a best fit peak at a lower 
redshift of z ~ 2.5. However, this is extremely sensitive to 
any evolution in the IR color distribution (but see § 5 to fol- 
low). We defer detailed fitting of our model to the available 
high redshift counts and FIRB to a companion paper (Lewis et 
al. 2003). Ultimately, redshift surveys of SIRTF and SCUBA 
galaxies will be required to test the detailed form of the evolu- 
tion function (see Chapman et al. 2003a for preliminary results 
on the SCUBA galaxy redshift distribution). 

5. COMPARISON WITH DATA OUT TO HIGHER REDSHIFTS 

Our comparison of evolving a single variable LF versus a bi- 
variate <^(C, C) (Figs. 6-8) clearly demonstrates that an over 
simplified representation of the galaxy population can miss 
large classes of detectable sources at high redshifts. It is then 
crucial to understand whether the IR color distribution itself has 
evolved with redshift. 

While much attention has been paid to the luminosity evo- 
lution of ULIRGs (e.g. Blain et al. 1999a,b; Chary & Elbaz 
2001; Malkan & Stecker 2001; Rowan-Robinson 2001; Chap- 
man et al. 2002b; Ivison et al. 2002), there has been no study 
of the TZ(60, 100) evolution. The simplest scenario, which we 
considered above, is that the 7^(60, 100) correlation with Ltir 
continues to higher luminosities, extrapolating the expressions 
derived in this work for the local TZ(60, 100) distribution. 

Existing data sets are largely inadequate for characterization 
of the high-z IR color distribution. ISOPHOT and SCUBA 
galaxies do not generally have accurate redshifts or even de- 
tections at sufficient wavelengths to measure the rest frame 
7^(60, 100). However, a survey of moderate redshift, luminous 
IRAS galaxies (Stanford et al. 2000) provides an initial foray 
into the higher redshift 7^(60, 100) distribution (for a comple- 
mentary analysis of this population, see also Blain, Barnard & 
Chapman 2003). This section will focus on the Stanford et al. 
IRAS galaxies, as well as the microjansky radio sources (Chap- 
man et al. 2003b) to assess our $(£, C) characterization. 



5.1. 60 ixm-selected LIRGs/ULIRGs at0.l<z< 0.9 

A recent sample of distant LIRGs and ULIRGs from Stan- 
ford et al. (2000) has an effective sensitivity only 2-3 times 
higher than that of the deepest ISOPHOT 170/im surveys, and 
we can use this sample to study the higher redshift <&(£, C). The 
Stanford et al. sample was selected from a positional cross- 
correlation of the IRAS Faint Source Catalog with the FIRST 
database. Objects from this set were selected for spectroscopy 
by virtue of following the well-known star-forming galaxy cor- 
relation between 1.4 GHz and 60 /im flux, and by being opti- 
cally faint on the POSS. Optical identification and spectroscopy 
were obtained for 108 targets at the Lick Observatory 3 m tele- 
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Fig. 9. — upper panel: The dependence of rest frame 7?.(60, 100) (roughly dust temperature) on Lfir is shown. LlRGs and ULIRGs from the Stanford et al. (2000) 
FlRST/60/^m selected sample, with redshifts from 0.1 to 0.9 (magenta circles), appear as a consistent extension to local IRAS-BGS galaxies (green circles) Many 
Stanford et al. sources are detected at 60 and 100 ^m (solid circles), although some have only tentative 100 /^m measurements when analysed in XSCANPI (open 
circles); those without 100 /^m detections at all are shown as limits. The fit to the Stanford et al. sources between 11. 5< logLFiR<12.5 is shown with lighter solid 
and dashed lines, for the ±lo" eiwelope. Two FIRBACK 170 /^m selected galaxies from Chapman et al. (2002a) are shown as large circles. Our best fit coirelation 
from section 3 and Icr deviations are overlaid as the thick solid and dashed lines. Arp220, the ULIRG typically used for comparison with high-z luminous galaxies, 
is identified with a star, lower panel: The 7?.(60, 100) ratio is shown as a function of redshift for the same sources. Overlapping luminosity regions can be compared 
directly. Triangles have logLpiRr^ll.SitO.S. Fits to the population overlapping in luminosity at Ltir=3x lO" Lq signalled by oversize triangles, is shown with a 
solid line, along with la, 2a distributions as dashed lines. The slopes of the fit line is not significantly different from zero. 



scope. Most objects show spectra typical of starburst galaxies, 
and do not show the high ionization lines of active galactic nu- 
clei. The redshift distribution covers 0.1 < z < 0.9, with 13 
objects at z > 0.5 and an average redshift of z = 0.3 1 . 

5.1.1. The dependence ofTZ(60, 100) on C 

We first address the issue of the dependence of width in C 
(= n(60, 100)) on the luminosity. Fig. 9a shows the 7e(60, 100) 
versus L distribution of the Stanford et al. sample relative to 



our fit C relation from section 3, overlaying the local BGS 
IRAS galaxies observed with SCUBA by Dunne et al. (2000) 
for reference. The rest frame TZ(6Q, 100) has been calculated 
by fitting the Dale et al. (2001, 2002) template SEDs to the ob- 
served frame 6Qfim, 100/^m, and radio points, and extracting 
the restframe 6Qfim, 100 /im fluxes. The moderate-z sources 
from Stanford et al. (2000) are not, however, all detected at 
100/im. We have run the IPAC application XSCANPI over all 
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103 Stanford sources to derive improved estimates of the 60 
and 100/im fluxes, finding 46 with spurious lOO/im measure- 
ments. The sources without 100 /xm detections are presented 
as lower Umits to 7?.(60, 100), while those with only marginal 
100/im detections are differentiated with open symbols. The 
fits to the distribution are then derived, using only sources with 
secure lOO/im detections. The fit result is shown as the lighter 
solid line with the ±1(7 lines shown as Ughter dashed Unes in 
the upper panel of Fig. 9. The sources with 100 /zm limits ap- 
pear to be uniformly distributed amongst the 100/zm detected 
sources, and a large systematic skew is unlikely. 

The median relation between Lfir and 7?.(60, 100) is similar 
to that observed locally, although it falls to the -Icr deviation of 
the local relation for the most luminous Stanford sources. How- 
ever, these most luminous sources are not very numerous and 
contribute negUgibly to the fit fines. Moreover, the fit did not 
take into account the lower limits to TZ(60, 100). The median 
of this luminous tail of Stanford sources on their own is in fact 
more consistent with the local relation. 

The width of the Stanford distribution matches the extrapo- 
lation from the 1.2 Jy envelope at Lfir=3.1 x 10^^ (Fig. 9a). 
However, the width of the Stanford sources exceed the local dis- 
tribution by a factor of ^50% at Lfir=4.6x 10^' Lq, and there 
is an apparent excess of LIRG/ULIRG class sources as cold as 
the two FIRBACK galaxies discussed in the introduction, FB 1- 
40 and FB 1-64 at their respective FIR luminosities (see §5.1.3). 
As the Stanford et al. sources were selected in the 60 /xm band, 
a direct comparison with our 1.2 Jy catalog is justified, with 
the assumption that the fainter flux limits of the Stanford et al. 
sample simply allows a probe of the higher luminosity distribu- 
tion out to moderate redshifts where evolution is still not likely 
to be a large effect. However, even for sources with detected 
100/im fluxes, the signal-to-noise ratio is not high, and photo- 
metric scatter may affect the width of the diagram. A further 
complication is that while the rest-frame 7^(60, 100) has been 
demonstrated to correlate well with cold dust temperature us- 
ing 850 /xm measurements (Dunne et al. 2000), the observed 
7^(60, 100) becomes less effective at constraining the cold tem- 
perature properties as redshifts increase beyond z ~ 0.5. The 
majority of the Stanford sources are at z ^ 0.4 and this should 
not result in a large effect on the average sample properties. 
Observations at 850/im will be required to properly assess the 
contribution of the Stanford sources at higher redshifts to the 
width of the distribution. 

We therefore conclude that the direct extrapolation from the 
local IRAS galaxies is consistent with the width of the Stanford 
et al. sources, with a possibility of an increasing width at lower 
luminosities in the Stanford sample. We note, however, that the 
width of even the local distribution at the highest luminosities 
is highly uncertain, and a change in the slope of the broadening 
function could be in agreement to within the errors. 

5.1.2. The dependence of 71(60, 100) on redshift 

We now address the possible evolution in the (>C,C) rela- 
tion with redshift. In Fig. 9b, the 7?,(60, 100) ratio is shown 
as a function of redshift for the same sources. The overlapping 
luminosity regions can be compared directly. Triangles have 
log Lfir/L0~1 1.5±0.5, in both IRAS-BGS and Stanford et al. 
galaxies. A fit to the overlapping luminosity populations (sig- 
nalled by oversize symbols) is shown with a solid line, along 
with Icr and 2a envelopes. 

The salient point of Fig. 9b is that the Stanford sources do not 
appear significantly different in the median than the local com- 



parison sample. The large fraction of Stanford sources with 
lower limits and uncertain 100/im fluxes may be responsible 
for the median point lying somewhat lower than local galaxies. 
The current data therefore do not support a scenario in which 
the median 7?.(60, 100) in C varies significantly with redshift out 
to z ^ 1 . The width in C of both luminosity bins is similar in 
Fig. 9b, and wider than the local distribution. This is consistent 
with the broadening of C with luminosity, and suggests that the 
redshift evolution affects the width in C only in so far as higher 
redshift sources have characteristically higher L. However, the 
statistics for sources with z > 0.5 are poor. With identified red- 
shifts for increased numbers of sub-mm galaxies (see Chapman 
et al. 2003a), this issue will be directly addressable. 

5.1.3. ISOPHOT 1 70 /xm galaxies at intermediate redshift 

Recently, Keck spectra and UKIRT high spatial resolu- 
tion near-IR imagery for two of the proposed highest redshift 
sources from the FIRBACK-Nl 170/xm survey were obtained 
(Chapman et al. 2002a). These authors found that the red- 
shifts of counterparts to the 170 /xm sources confirm that both 
sources are ULIRGs, but that their redshifts are significantly 
lower than implied by fitting a typical ULIRG SED to their 
farlR/sub-mm/radio SEDs. This indicates that they have cooler 
dust temperatures, T^ ^ 30 K, than the canonical ULIRG values 
(Trf ~ 45 K). As previous models have failed to predict the ex- 
istence of such cold and luminous sources, these two galaxies 
therefore confirm the importance of considering the full $(£, 
C) in evolutionary models. 

In Fig. 9, we overplot these 2 FIRBACK ULIRGs on the lo- 
cal IRAS BGS and the Stanford LIRG/ULIRG sample; FNl-40 
and FNl-64 are cold objects, lying near the -2a level below 
the local median relation. While FBI -40,64 appear as cold ob- 
jects relative to the local distribution in $(£, C), there is not 
yet enough followup to ISOPHOT sources to know if C has 
substantially evolved out to z ~ 1. FN 1-40,64 remain at the pe- 
riphery of the distribution, but are certainly not such unusual 
objects when compared to similar redshift ULIRGs from the 
Stanford et al. sample. 

5.2. The C-C relation for microjansky radio sources 

In this subsection, we test the lower luminosity C-C relation 
out to z ^ 0.6 using a sample of radio galaxies with 20cm fluxes 
ranging from 30/LtJy to 500/iJy from Chapman et al. (2003b). 
These radio sources were observed in the sub-mm, and all have 
spectroscopic redshifts ranging from z=0. 1-0.6. From this in- 
formation we can estimate Td and Lfir by assuming the FIR- 
radio correlation for star forming galaxies (Helou et al. 1985) 
and a synchrotron spectrum with an index of a = -0.8 (Richards 
2000). T(] is estimated by taking the SED template from Dale 
& Helou (2002) which has the appropriate temperature to fit 
the radio and sub-mm points at the fixed redshift of the source. 
Note that we are only indirectly testing the C-C relation with 
these galaxies (we are estimating both C and C through local 
empirical correlations), and we cannot deconvolve viiriations in 
the C-C relation from the FIR-radio correlation. 

We then restrict the sample to those sources with Lfir< 
10^' L0 as the spectroscopic incompleteness introduces a se- 
vere bias for more luminous (and typically optically fainter) 
galaxies (as discussed in Chapman et al. 2003b). This also re- 
duces the radio-bright AGN contribution to the sample, which 
will have much larger apparent Lfir calculated in this manner 
due to the AGN generated radio excess. 
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Fig. 10. — Dust temperature (Tj) versus the log of FIR luminosity for ra- 
dio sources with Lfir< 10"Lq and spectroscopic redshifts. Sources span 
z = 0.1-0.6, and have been divided into three equal number bins. The Tj 
and FIR values have been calculated from the redshift and the sub-mm/radio 
measurements as described in the text. Our derivation of the median and in- 
terquartile range of local IRAS galaxies from the 1 .2 Jy catalog are shown. The 
agreement of the higher redshift sources with the local relation is remarkable. 



In Fig. 10 we plot Td versus FIR luminosity, along with our 
derived local IRAS galaxy relation from § 3. The average points 
in Fig. 10 do not include sources individually detected in the 
submillimeter, and the variance from measurement errors is ex- 
pected to be comparable to that from scatter in the dust proper- 
ties. Thus our points can only be used as a characterization of 
the average C-C relation. 

Fig. 10 demonstrates that the Lfir< IO^'L© radio sources, 
spanning the redshift range z = 0. 1-0.6 appear to follow the lo- 
cal IRAS color distribution very well. This suggests jointly that 
neither the C-C relation nor the FIR-radio correlation can de- 
viate significantly from local values for these sources. The 
analysis provides a direct consistency check of the lower lumi- 
nosity tail of local IRAS galaxies with higher redshift star form- 
ing galaxies, of which the microjansky radio population should 
represent higher redshift specimens (e.g., Richards 2000). 



CONCLUSIONS 

We have analyzed the 7^(60, 100) distribution of local IR- 
luminous galaxies finding a best fit, low order analytical ex- 
pression for the bi-variate luminosity function, C). A 
log-normal distribution about a dual power law in the median 
7^(60, 100) versus Ltir is demonstrated to provide the best fit 
of the IRAS population over 4 orders of magnitude in Ljir. 

We then studied the redshift evolution of <!>(£, C), using a 
luminosity evolution paradigm. We demonstrated that while 
similar luminosity regimes are detectable at 24 /im and 850 /im, 
the dust temperatures represented by the detectable sources can 
differ by factors of > 2. For a flux limited survey, the bi-variate 
$(£, C) model predicts about twice as many warm galaxy de- 
tections as the equivalent single variable model, even though 
the underlying luminosity distributions are the same in the two 
models. In a similar manner, recently discovered populations 
of cold, luminous galaxies can be predicted naturally within our 
model. Consideration of the model differences also addresses 
the question of the SCUBA/SIRTF galaxy overlap, a strong 
function of the dominant dust temperature. A cold (Tj < 30 K) 
SCUBA population may be difficult to detect even in the deep- 
est SIRTF exposures, as revealed in the high redshift excess of 
cold SCUBA sources compared with the SIRTF population. 

We compare our derived color relation with existing data for 
higher redshift IR galaxies. General agreement with the local 
color relation is found amongst all higher redshift galaxy pop- 
ulations considered. We therefore find no significant evidence 
for a variation in the median 7^(60, 100) ratio as a function of 
redshift. 

We do find tentative evidence for a broadening of the IR color 
distribution at higher C from both our local 1 .2 Jy sample, and 
the moderate redshift high luminosity sample of Stanford et al. 
(2000). This suggests that redshift evolution does not in itself 
affect the width in C. The baseline assumption that the locally 
measured C-C relation holds over all redshifts is therefore a 
resonable working hypothesis, to be verified with SIRTF data. 
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